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With the goal of developing predictive ab-initio capability for light and medium-mass nuclei, two- 
nucleon and three-nucleon forces from chiral effective field theory are optimized simultaneously to 
low-energy nucleon-nucleon scattering data, as well as binding energies and radii of few-nucleon 
systems and selected isotopes of carbon and oxygen. Coupled-cluster calculations based on this 
interaction, named NNLOsat, yield accurate binding energies and radii of nuclei up to '^'^Ca, and 
are consistent with the empirical saturation point of symmetric nuclear matter. In addition, the 
low-lying collective — 3“ states in and ^°Ca are described accurately, while spectra for 
selected p- and sd-shell nuclei are in reasonable agreement with experiment. 

PACS numbers: 21.30.-x, 21.10.-k, 21.45.-v, 21.60.De 


Introduction - Interactions from chiral effective field 
theory (EFT) [IHl] and modern applications of renor¬ 
malization group techniques have opened the door 
for a description of atomic nuclei consistent with the un¬ 
derlying symmetries of quantum chromodynamics, the 
theory of the strong interaction. Ghiral nuclear forces 
can be constructed systematically from long-range pion 
physics augmented by contact interactions. Over the past 
decade, the renaissance of nuclear theory based on real¬ 
istic nuclear forces and powerful computational methods 
has pushed the frontier of ab initio calculations from few- 
body systems and light nuclei HI HI HOl to medium-mass 
nuclei 

One of the main challenges in ab-initio calculations 
is the accurate [20j reproduction of binding energies and 
radii of finite nuclei simultaneously with the empirical nu¬ 
clear matter saturation point (binding energy per nucleon 
E/A « 16MeV at Fermi momentum kp ~ 1.33fm“^) 
and incompressibility (250 < Kq < 315MeV [H]). 
For instance, lattice EFT calculations at next-to-next-to 
leading order (NNLO) employ a phenomenological four- 
nucleon contact force (a correction beyond NNLO) to 
counter the overbinding in nuclei heavier than ^^G HZ], 
while the radii of and are still too small [^I23| . 
Ab initio calculations overbind medium-mass and heavy 
nuclei by about 1 MeV per nucleon, underestimate charge 
radii [Hj, and yield too large separation energies [ZSj. 
The status of chiral-force predictions for binding ener¬ 
gies and charge radii in finite nuclei is summarized in 
Fig. a with dark grey symbols representing the predic¬ 
tions of various state-of-the-art calculations. This is a 
serious shortcoming of current chiral Hamiltonians as it 
prevents theory from making accurate predictions when 


extrapolating to higher masses. The problem with the 
reproduction of nuclear matter saturation properties has 
been discussed extensively in the literature and 

various solutions have been proposed, ranging from short- 
range correlations and Pauli blocking effects to the inclu¬ 
sion of many-body forces. 
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FIG. 1. (Color online) Ground-state energy (negative of bind¬ 
ing energy) per nucleon (top), and residuals (differences be¬ 
tween computed and experimental values) of charge radii 
(bottom) for selected nuclei computed with chiral interac¬ 
tions. In most cases, theory predicts too small radii and too 
large binding energies. References: a |40l I41| . b [^, c |23| . 
d [ 22 ], e [42], / [ 43 ], g [44], h [45], i [46]. The red diamonds 
are NNLOsat results obtained in this work. 

We start from the optimization of the chiral interac¬ 
tion at NNLO. Traditionally, one takes the pion-nucleon 
coupling constants c/s either from pion-nucleon scat- 
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tering im or from peripheral partial waves in the 
nucleon-nucleon (NN) sector 150] , while the remain¬ 
ing coupling constants (denoted as low-energy constants 
(LECs)) are adjusted in the NN sector. The correspond¬ 
ing optimizations consider scattering data up to lab¬ 
oratory energies of Tuab ~ 350 MeV. In a subsequent 
step, the remaining LECs of the leading three-nucleon 
(NNN) forces [5TH5^ are adjusted to data on A < 4 
systems For details, we refer the reader to 

Refs. Hitherto such a strategy has not produced 

interactions that simultaneously describe bulk properties 
of both nuclei and nuclear matter |56| . 

Our optimization strategy is based on a different ap¬ 
proach. Most importantly, we optimize NNN forces 
together with NN forces. This is consistent with the 
idea of EFT that improvements are made order by or¬ 
der and not nucleon by nucleon. The simultaneous opti¬ 
mization of NN and NNN forces is important because 
the long-range contributions of the NNN force contain 
LECs from pion-nucleon vertices that also enter the NN 
force. Moreover, in addition to low-energy NN data and 
the binding energies and charge radii of ^H, our 

set of fit-observables also contains data on heavier nuclei; 
namely, binding energies and radii of carbon and oxygen 
isotopes. This is a major departure from the traditional 
approach that seeks to adjust LECs to data on few-body 
systems with H = 2, 3,4 and then attempts to extrapo¬ 
late to nuclei with A ^ 1 and to infinite nuclear mat¬ 
ter. The following arguments motivate the strategy of 
including heavier nuclei into the optimization: First, no 
reliable experimental data constrain the isospin T = 3/2 
components of the NNN force in nuclei with mass num¬ 
bers H = 3,4 (see Refs. [571 HH] for more discussion and 
prospects). Second, since our goal is to describe nu¬ 
clear properties at low energies, LECs are adjusted to 
low-energy observables (as opposed to the traditional ad¬ 
justment of two-nucleon forces to NN scattering data at 
higher energies). Third, the impact of many-body effects 
entering at higher orders (e.g., higher-rank forces) might 
be reduced if heavier systems, in which those effects are 
stronger, are included in the optimization. 

Besides these theoretical arguments, there is also one 
practical reason for a paradigm shift: predictive power 
and large extrapolations do not go together. In tradi¬ 
tional approaches, where interactions are optimized for 
H = 2,3,4, small uncertainties in few-body systems (e.g., 
by forcing a rather precise reproduction of the H = 2, 3,4 
sectors at a rather low order in the chiral power count¬ 
ing) get magnified tremendously in heavy nuclei, see for 
example Ref. [23] . Consequently, when aiming at reliable 
predictions for heavy nuclei, it is advisable to use a model 
that performs well for light- and medium-mass systems. 
In our approach, light nuclei are reached by interpola¬ 
tion while medium-mass nuclei by a modest extrapola¬ 
tion. In this context, it is worth noting that the most 
accurate calculations for light nuclei with A < 12 [51] 


employ NNN forces adjusted to 17 states in nuclei with 
H < 8 [60]. Finally, we point out that nuclear satu¬ 
ration can be viewed as an emergent phenomenon. In¬ 
deed, little in the chiral EFT of nuclear forces suggest 
that nuclei are self-bound systems with a central den¬ 
sity (or Fermi momentum) that is practically indepen¬ 
dent of mass number. This viewpoint makes it prudent 
to include the emergent momentum scale into the opti¬ 
mization, which is done in our case by the inclusion of 
charge radii for ^H, ^djje, and ^®0. This is similar 
in spirit to nuclear mean-field calculations m and nu¬ 
clear density functional theory [621163] where masses and 
radii provide key constraints on the parameters of the 
employed models. 

Optimization protocol and model details - We seek to 
minimize an objective function to determine the opti¬ 
mal set of coupling constants of the chiral NN-\-NNN 
interaction at NNLO. Our dataset of fit-observables in¬ 
cludes the binding energies and charge radii of ^H, 

and as well as binding energies of 22,24,25 q g^g 
summarized in Table [T] To obtain charge radii rch from 
computed poin-proton radii Cpp we use the standard ex¬ 
pression [Mj: = {rip) + (i?p) + f (i?n) + 

where = 0.033 fm^ (Darwin-Foldy correction), 

Rl = -0.1149(27) fm2 [ 55 ] and Rp = 0.8775(51) fm [ 55 ]. 
In this work we ignore the spin-orbit contribution to 
charge radii m- From the NN sector, the objec¬ 
tive function includes proton-proton and neutron-proton 
scattering observables from the SM99 database [ 55 ] up 
to 35 MeV scattering energy in the laboratory system as 
well as effective range parameters, and deuteron proper¬ 
ties (see Table 0- The maximum scattering energy was 
chosen such that an acceptable fit to both NN scattering 
data and many-body observables could be achieved. 


TABLE I. Binding energies (in MeV) and charge radii (in 
fm) for ^H, ^’'*He, and i6>22,23,24,25q gjgpjQygfj jg ^pg gp. 

timization of NNLOaat- 



Ega 

Exp. [69] 

^ch 

Exp. [65] [66] 


8.52 

8.482 

1.78 

1.7591(363) 

®He 

7.76 

7.718 

1.99 

1.9661(30) 

"He 

28.43 

28.296 

1.70 

1.6755(28) 

14c 

103.6 

105.285 

2.48 

2.5025(87) 

I 60 

124.4 

127.619 

2.71 

2.6991(52) 

22 o 

160.8 

162.028(57) 



24o 

168.1 

168.96(12) 



250 

167.4 

168.18(10) 




In the present optimization protocol, the NNLO chiral 
force is tuned to low-energy observables. The comparison 
with the high-precision chiral NN interaction N^LOem 
[ 45] and experimental data presented in Table |n] demon¬ 
strates the quality of NNLOsat at low energies. 

The results for and ^djjg (and ®Li) were computed 
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TABLE II. Low-energy NN data included in the optimiza¬ 
tion. The scattering lengths a and effective ranges r are in 
units of fm. The proton-proton observables with superscript 
C include the Coulomb force. The deuteron binding energy 
{Ed, in MeV), structure radius {rn, in fm), and quadrupole 
moment {Qd, in fm^) are calculated without meson-exchange 
currents or relativistic corrections. The computed d-state 
probability of the deuteron is 3.46%. 



NNLOsat 

N^LOem [ 49 ] 

Exp. 

Ref. 

CLpp 

-7.8258 

-7.8188 

-7.8196(26) 

|70| 

1 pp 

2.855 

2.795 

2.790(14) 

EOj 

^nn 

-18.929 

-18.900 

-18.9(4) 

[71] 

^ Tin 

2.911 

2.838 

2.75(11) 

m 

f^np 

-23.728 

-23.732 

-23.740(20) 

m 

T'np 

2.798 

2.725 

2.77(5) 

m 

Ed 

2.22457 

2.22458 

2.224566 

m 

rD 

1.978 

1.975 

1.97535(85) [71] 

Qd 

0.270 

0.275 

0.2859(3) 

[73] 


with the no-core shell model (NCSM) 0 [TU] accompa¬ 
nied by infrared extrapolations [75]. The NNN force 
of NNLOgat yields about 2 MeV of binding energy for 
“^He. Heavier nuclei are computed with the coupled- 
cluster method (see Ref. [75] and the discussion below). 

A total of 16 LECs determine the strengths of the NN 
contact potential, the ttN potential in the NN+NNN 
sector, and the NNN contacts. The LECs are con¬ 
strained simultaneously by the optimization algorithm 
POUNDerS [B5]. We employ standard nonlocal reg¬ 
ulators in the construction of the potential, see e.g.. 
Refs 112] for details. This type of regulator improves 
the convergence of nuclear matter calculations |56j . In 
detail, the regulator functions consist of exponentiated 
Jacobi momenta p divided by a cutoff value A, i.e., 
~ exp[(p/A)^”]. For the present work, we set n = 3 
and A = 450 MeV. Furthermore, the subleading two- 
pion exchange in the NN interaction is regularized us¬ 
ing spectral function regularization with a cutoff Asfr = 
700 MeV. The details of this procedure can be found in 
Refs. [711 US] . 

The objective function is numerically expensive, re¬ 
quiring us to adopt some approximations when comput¬ 
ing nuclei with A > 4. In the optimization, we employed 
a model space of 9 oscillator shells for the NN interac¬ 
tion, the energy cutoff Esmax = 8hV, for the NNN forces, 
and the coupled-cluster method in its singles and dou¬ 
bles approximation (CCSD). We use nucleus-dependent 
estimates for larger model spaces and triples-cluster cor¬ 
rections based on Ref. m- During the optimization, we 
verified that these estimates were accurate by perform¬ 
ing converged calculations. In our final computation of 
the objective function and for the results presented in this 
paper, we employ much larger model spaces and coupled- 
cluster methods with higher precision. 

The coupled-cluster calculations are based on the in¬ 


trinsic Hamiltonian H = T — Tcm + Vnn + R/v wtv to min¬ 
imize spurious center-of-mass effects [15] [791 [80] . For the 
binding energies presented in this paper we employ the 
A-CCSD(T) approximation [TSIIHIIIHI] in a model space 
consisting of 15 oscillator shells with HU = 22 MeV. The 
NNN forces are limited to an energy cutoff Esmax = 
16hn, and truncated at the normal-ordered two-body 
level in the Hartree-Fock basis [SS] |84] . We also include 
the leading-order residual NNN contribution to the total 
energy as a second-order perturbative energy correction 
[56] ■ computed with Famax = 

To compute excited states in, and around, nuclei with 
closed shells, we employ equation-of-motion coupled- 
cluster methods [MllSHUS^ : these are accurate for excited 
states that are generalized particle-hole excitations of low 
rank. For instance, is computed with the charge- 
symmetry breaking equation-of-motion method from the 
closed sub-shell nucleus see Ekstrdmetal. [89]. Sim¬ 
ilar comments apply to 22,24p^ ppg intrinsic charge radii 
are computed from the two-body density matrix (2BDM) 
in the CCSD approximation [9D]. Benchmark calcula¬ 
tions of the '^He charge radius shows that the 2BDM re¬ 
sult is 1% larger than the NCSM result. Intrinsic charge 
densities are computed using the one-body density ma¬ 
trix and correcting for the Gaussian center-of-mass wave- 
function [Tiini]. In the case of this approach has 
been validated against 2BDM to four significant digits. 

The values for the LECs that result from the optimiza¬ 
tion and define the chiral potential NNLOsat are listed in 
Table [hi] We note that the pion-nucleon LECs ci, C 3 and 


TABLE III. The values of the LECs for the NNLOsat inter¬ 
action. The constants d, Ci, and Ci are in units of GeV“^, 
lO^GeV”^, and lO'^GeV”'*, respectively. 


LEG 

Value 

LEG Value 

LEG Value 

Cl 

-1.12152120 

C3 -3.92500586 

C4 3.76568716 


-0.15814938 Cr? -0.15982245 
*^0 

Cig -0.15915027 


2.53936779 

C^s^ 1.00289267 

C3g^ -0.17767436 

Cl Pi 

0.55595877 

Csp^ 1.39836559 

Csp^ -1.13609526 

Cssj-Soj 

Ce 

0.60071605 

-0.03957471 

Csp^ -0.80230030 

CD 0.81680589 


C 4 are in the range of the published values [IS] HU | 22 ] ■ 
Following Ref. [12], we set the pion-decay constant = 
92.4 MeV and the axial-vector coupling constant gA = 
1.29. The value for gA is greater than the experimental 
estimate gA = 1.276 [93] to account for the Goldberger- 
Treiman discrepancy. We use the following neutron, pro¬ 
ton, and nucleon masses: = 939.5653 MeV, rrip = 

938.272 MeV, and rriN = 938.9184 MeV, respectively. 
For the pion masses we used m^± = 139.5702 MeV and 
m^o = 134.9766 MeV. For the NN scattering data up to 
35 MeV a total y^/datum « 4.3 was reached. Represen¬ 
tative phase shifts are shown in Fig.[^ The phase shifts 
at higher scattering energies, demonstrates that NNLOsat 
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is at the limits of expectations one can have for an inter¬ 
action at this chiral order. Furthermore, the accuracy of 
NNLOsat in the few-body sector is similar to other chiral 
interactions at order NNLO [511 IHl] ■ 
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FIG. 2. (Color online) Selected neutron-proton scattering 
phase-shifts as a function of the laboratory scattering energy 
Thab- (Top) NNLOsat prediction (solid lines) compared to 
the Nijmegen phase shift analysis [95] (symbols) at low ener¬ 
gies Tbab < 35 MeV. Note the two vertical scales. (Bottom) 
Neutron-proton scattering phase shifts from NNLOsat (red di¬ 
amonds) compared to the Nijmegen phase shift analysis (black 
squares) and the NNLO potentials (green) from Ref. | 77| . 


Predictions — We begin with predictions for the /?- 
decay half-life of ^H. The reduced matrix element 
|(^He||i?j^|pH)| = 0.6343 compares well to the corre¬ 
sponding experimental value of 0.6848 ± 0.0011 [551196] . 
Figure [l] shows that binding energies and charge radii of 
the p-shell nuclei ®He, ^'^C, and are in good agree¬ 
ment with experiment. For ®He the computed binding 
energy and charge radius are 30.9 MeV and 1.91 fm, re¬ 
spectively, and in good agreement with the experimental 
binding energy 31.5 MeV |69j and experimental charge 
radius 1.959(16) fm [57]. For ®’®Li we compute a binding 
energy of 32.4(4) MeV and 43.9 MeV, respectively, which 
compare well with experiment (32.0 MeV and 45.34 MeV 


[69]). The charge radius of ®Li with NNLOsat is 2.22 fm, 
also consistent with the measured value of 2.217(35) fm 

[98] . We now discuss results for excited states in ®Li, 
14c, ^^N, and see Fig.[^ The nucleus ®Li is difficult 
to compute because it is bound by only 1.5 MeV relative 
to the threshold for deuteron emission. Effects of contin¬ 
uum is expected to lower the 2+ resonances significantly 

[99] . thus we conclude that our results are in reason¬ 
able agreement with experiment. We also compared the 
spectra computed in the NCSM and agreement with the 
coupled-cluster prediction is good. The binding energy 
computed from two-particle attached equation-of-motion 
method [80l|87] is 30.9MeV and in reasonable agreement 
with the NCSM extrapolated result 32.4(4) MeV. Our 
predictions for the excited states of and ^'^N agree 
with experiment except for the ij state in ^^N. 



FIG. 3. (Color online) Energies (in MeV) of selected excited 
states for various nuclei using NNLOsat. For ®Li we also 
include spectra from the NCSM (dotted lines), and isospin 
quantum numbers are also given. The NCSM results were 
obtained with Nmax = 10 and fill = 16 MeV. Parenthesis de¬ 
note tentative spins assignments for experimental levels. Data 
are from Refs. [TOMnSj . 


The ab-initio computation of negative-parity states in 
particularly the 3)" state at 6.13 MeV [23 [Ml [Ml 
[TT)5] has been a long-standing theoretical challenge. We 
computed this state, dominated by about 90% of Ip-lh 
{Pi /2 '^ 5 / 2 ) excitations, at 6.34 MeV. The energy of the 

3j" state is strongly correlated with the charge radius of 
with smaller charge radii leading to higher excitation 
energies. For Ip-lh excited states, the excitation energy 
depends on the particle-hole gap and therefore on one- 
nucleon separation energies of the A = 16 and A = 17 
systems. The charge radius depends also on the proton 
separation energy Sp. For we find Sp = 10.69MeV 
and the neutron separation energy S'„(^’^0) = 4.0 MeV, 
in an acceptable agreement with the experimental values 
of 12.12 MeV and 4.14 MeV, respectively. For ^^F we find 
Sp = 0.5 MeV, to be compared with the experimental 
threshold at 0.6 MeV. 
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The inset of Fig. shows that the 2^ state in 
also comes out well, suggesting a \p-lh nature. However, 
the state is about 1.5 MeV too high compared to 
experiment. This state is dominated by Ip-lh excitations 
from the occupied pi /2 to the unoccupied Si /2 orbitals. In 
the 1/2+ state is computed at an excitation energy of 
2.2 MeV, which is about 1.4 MeV too high. This probably 
explains the discrepancy observed for the 1“ state in 
Figure shows that the experimental charge-density 
of is well reproduced with NNLOsat, and our charge 
form factor is, for momenta up to the second diffrac¬ 
tion maximum, similar in quality to what Mihaila and 
Heisenberg m achieved with the Avl8 -|- UIX poten¬ 
tial. For the heavier isotopes 22,24 q 22,24p pjg 

shows good agreement between theory and experiment 
for excited states. For ^^F our computed spin assign¬ 
ments agree with results from shell-model Hamiltoni¬ 
ans |106j and with recent ab initio results [89] . The bind¬ 
ing energies for 22.24p 103.7 MeV, 163 MeV and 

175.1 MeV, respectively, in good agreement with data 
(104.7 MeV, 167.7 MeV and 179.1 MeV). We also com¬ 
puted the intrinsic charge (matter) radii of 22.24 q 
obtained 2.72 fm (2.80 fm) and 2.76 fm (2.95 fm), re¬ 
spectively. The matter radius of agrees with the 
experimental result from Kanungo et al. m- We note 
that the computed spectra in is too compressed com¬ 
pared to experiment (theory yields 0.7 MeV compared 
to 1.9 MeV for the first excited 2+ state), possibly due 
to the too high 1/2+ excited state in In general, 
the quality of our spectra for sd-shell nuclei is compara¬ 
ble to those of recent state-of-the-art calculations with 
chiral Hamiltonians jUJ I1()7H1()9| , while radii are much 
improved. 



FIG. 4. (Color online) Charge density in computed as in 
Ref. [llOj compared to the experimental charge density [lllj . 
The inset compares computed low-lying negative-parity states 
with experiment. 

For ^°Ca the computed binding energy E = 326 MeV, 
charge radius Teh = 3.48 fm, and £'(3/) = 3.81 MeV 
all agree well with the experimental values of 342 MeV, 


3.4776(19) fm |5S|, and 3.736 MeV respectively. We 
checked that our energies for the 3/ states in and 
^'^Ca are practically free from spurious center-of-mass ef¬ 
fects. The results for "‘^Ca illustrate the predictive power 
of NNLOsat when extrapolating to medium-mass nuclei. 

Finally, we present predictions for infinite nuclear mat¬ 
ter. The accurate reproduction of the saturation point 
and incompressibility of symmetric nuclear matter has 
been a challenge for ab initio approaches, with represen¬ 
tative results from chiral interactions shown in Fig. 
The solid line shows the equation of state for NNLOsat- 
Its saturation point is close to the empirical point, and its 
incompressibility K = 253 lies within the accepted em¬ 
pirical range [2T|. At saturation density, coupled-cluster 
with doubles yields about 6 MeV per particle in correla¬ 
tion energy, while triples corrections (and residual NNN 
forces beyond the normal-ordered two-body approxima¬ 
tion) yield another 1.5 MeV. 



kj (fm 

FIG. 5. (Color online) Equation of state for symmetric nu¬ 
clear matter from chiral interactions. Solid red line: pre¬ 
diction of NNLOsat- Blue dashed-dotted and black dashed 
lines: Ref. |56| . Symbols (red diamond, blue circle, black 
square) mark the corresponding saturation points. Trian¬ 
gles: saturation points from other models (upward triangles 
[33], rightward triangles dm, downward triangles [36] )■ The 
corresponding incompressibilities (in MeV) are indicated by 
numbers. Green box: empirical saturation point. 

Let us briefly discuss the saturation mechanism. Sim¬ 
ilar to Vowfe potentials [S], the NN interaction of 
NNLOsat is soft and yields nuclei with too large bind¬ 
ing energies and too small radii. The NNN interactions 
of NNLOsat are essential to arrive at physical nuclei, sim¬ 
ilarly to the role oi NNN forces in the saturation of nu¬ 
clear matter with low-momentum potentials [33] . This 
situation is reminiscent of the role the three-body terms 
play in nuclear density functional theory [1 13] . 

Summary - We have developed a consistently opti¬ 
mized interaction from chiral EFT at NNLO that can be 
applied to nuclei and infinite nuclear matter. Our guide¬ 
line has been the simultaneous optimization of W and 
NNN forces to experimental data, including two-body 
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and few-body data, as well as properties of selected light 
nuclei such as carbon and oxygen isotopes. The optimiza¬ 
tion is based on low-energy observables including binding 
energies and radii. The predictions made with the new 
interaction NNLOsat include accurate charge radii and 
binding energies. Spectra for “^^Ca and selected isotopes 
of lithium, nitrogen, oxygen and fluorine isotopes are well 
reproduced, as well as the energies of 3j" excitations in 
and ^°Ca. To our knowledge, NNLOgat is currently 
the only microscopically-founded interaction that allows 
for a good description of nuclei (including their masses 
and radii) in a wide mass-range from few-body systems 
to medium-mass. 
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